Section: New Results

Inverse problems for Poisson-Laplace equations

Participants : Laurent Baratchart, Sylvain Chevillard, Juliette Leblond, Konstantinos Mavreas, Christos Papageorgakis, Dmitry Ponomarev.

This section is concerned with inverse problems for 3-D Poisson-Laplace equations, among which source recovery issues. Though the geometrical settings differ in Sections 7.1.1 and 7.1.2 , the characterization of silent sources (those giving rise to a vanishing field) is one common problem to both which has been resolved in the magnetization setup [33] .

Inverse problems in medical imaging

This work is conducted in collaboration with Jean-Paul Marmorat and Nicolas Schnitzler, together with Maureen Clerc and Théo Papadopoulo from the Athena EPI.

In 3-D, functional or clinical active regions in the cortex are often modeled by pointwise sources that have to be localized from measurements taken by electrodes on the scalp of an electrical potential satisfying a Laplace equation (EEG, electroencephalography). In the works [38] [5] on the behavior of poles in best rational approximants of fixed degree to functions with branch points, it was shown how to proceed via best rational approximation on a sequence of 2-D disks cut along the inner sphere, for the case where there are finitely many sources (see Section 4.2 ).

In this connection, a dedicated software FindSources3D (see Section 3.4.2 ) is being developed, in collaboration with the team Athena and the CMA. We continued this year algorithmic developments, prompted by a fruitful collaboration with the firm BESA, namely automatic detection of the number of sources (which was left to the user until recently). It appears that, in the rational approximation step, multiple poles possess a nice behavior with respect to branched singularities. This is due to the very physical assumptions on the model (for EEG data, one should consider triple poles). Though numerically observed in [7] , there is no mathematical justification so far why multiple poles generate such strong accumulation of the poles of the approximants. This intriguing property, however, is definitely helping source recovery. It is used in order to automatically estimate the “most plausible” number of sources (numerically: up to 3, at the moment). Further, a modular and ergonomic platform version of the software is under development.

In connection with these and other brain exploration modalities like electrical impedance tomography (EIT), we are now studying conductivity estimation problems. This is the topic of the PhD research work of C. Papageorgakis (co-advised with the Athena project-team and BESA GmbH). In layered models, it concerns the estimation of the conductivity of the skull (intermediate layer). Indeed, the skull was assumed until now to have a given isotropic constant conductivity, whose value can differ from one individual to another. A preliminary issue in this direction is: can we uniquely recover and estimate a single-valued skull conductivity from one EEG recording? This has been established in the spherical setting when the sources are known, see [17] . Situations where sources are only partially known and the geometry is more realistic than a sphere are currently under study. When the sources are unknown, we should look for more data (additional clinical and/or functional EEG, EIT, ...) that could be incorporated in order to recover both the sources locations and the skull conductivity. Furthermore, while the skull essentially consists of hard bone part that may be assumed to have constant electrical conductivity, it also contains spongy bone compartments. These two distinct parts of the skull possess quite different conductivities. The influence of that second value on the overall model is now being studied [19] .

Inverse magnetization issues in the thin-plate framework

This work is carried out in the framework of the “équipe associée Inria” IMPINGE, comprising Eduardo Andrade Lima and Benjamin Weiss from the Earth Sciences department at MIT (Boston, USA) and Douglas Hardin, Michael Northington and Edward Saff from the Mathematics department at Vanderbilt University (Nashville, USA).

Localizing magnetic sources from measurements of the magnetic field away from the support of the magnetization has been the fundamental issue under investigation by IMPINGE. The goal was to determine magnetic properties of rock samples (e.g. meteorites or stalactites) from fine field measurements close to the sample that can nowadays be obtained using SQUIDs (superconducting quantum interference devices). Currently, rock samples are cut into thin slabs and the magnetization distribution is considered to lie in a plane, which makes for a somewhat less indeterminate framework than EEG as regards inverse problems because “less” magnetizations can produce the same field (for the slab has no inner volume). Note however that EEG data consist of values of the normal current and of the associated potential, while in the present setting only values of the normal magnetic field are measured.

Figure 3 presents a schematic view of the experimental setup: the sample lie on a horizontal plane at height 0 and its support is included in a rectangle. The vertical component $B_3$ of the field produced by the sample is measured on points of a horizontal $N\times N$ rectangular grid at height $h$.

Over the previous years, we mainly focused on developing techniques to recover magnetizations with rather sparse support. To this end, we set up a heuristic procedure to recover sparse magnetizations, based on iterative truncation of the support of the recovered magnetization. In this heuristics, magnetizations were represented by dipoles placed at the points of a regular rectangular $n\times n$, which seemed general enough a model class to correctly approximate the magnetizations commonly encountered in samples.

The procedure turned out to be poor when trying to recover the magnetization itself, due to the severe ill-posedness of the problem and the unexpected existence of magnetizations that produce almost no field at the height where measurements are performed, although the corresponding magnetic distributions strongly differ from truly silent distributions. Nevertheless, whenever the support could be significantly shrunk while keeping the error small (i.e., explaining the data satisfactorily), estimates of the net moment so far, based on the dipolar model obtained by inversion, have been good.

This suggests that recovering the net moment and recovering the magnetization are rather different problems, the first one being less ill-posed than the second. Although the information provided by the net moment of the sample seems to be much weaker than knowing the full magnetic distribution, its importance has been emphasized by the geophysicists at MIT for at least three reasons:

• It yields important geological information on the sample in particular to estimate the magnitude of the ambient magnetic field at the time the rock was formed.

• It can be estimated independently to some extent, using a magnetometer, thereby allowing one to cross-validate the approach.

• From a computation point of view, knowledge of the net moment should lead to numerically stable reconstruction of an equivalent unidirectional magnetization. The support of the latter would provide us with valuable information to test for unidirectionality of the true magnetization, which is an important question to physicists in connection with rocks history and formation.

This year, we addressed the problem of directly recovering the net moment, without recourse to full inversion. Indeed, the latter is rather inefficient as it requires using a cluster and even then, for some samples, days of evaluation in order to obtain only a coarse estimate of the net moment. This research effort led us to investigate three different and complementary approaches.

First, we improved over Fourier based techniques previously designed by reformulating the problem with the help of the Kelvin transform. This gave us an asymptotic expansion of the net moment involving, at the first order, the integrals $\int \int B_3(x,y,h)\mathrm{d}x\mathrm{d}y$, $\int \int x{B}_3(x,y,h)\mathrm{d}x\mathrm{d}y$ and $\int \int y{B}_3(x,y,h)\mathrm{d}x\mathrm{d}y$, computed on a disc with large radius. Although the method is promising, the computations are quite involved and we did not manage yet to obtain higher-order terms. This is a part of D. Ponomarev PhD work.

In parallel, and based on the results obtained with Fourier transform, we investigated a second approach, consisting in directly computing asymptotic expansions of the above integrals, on several domains (namely, the 2-D balls of radius $R$ for the 1, 2 and $\infty$ norm, that are squares, disks, diamonds). In all cases, we get

where $\langle m_1\rangle$ is the moment of the first component $m_1$ of the magnetization and $\langle t_1{m}_3\rangle$ is the first moment of $m_3$ with respect to the first variable. The constants $\alpha$ and $\beta$ depend on the domain where the integral is computed. Therefore, an appropriate linear combination of the integrals computed on the different domains allows us to compute $\langle m_1\rangle$ with an accuracy of $𝒪(1/R^3)$. Similar results are obtained for $\langle m_2\rangle$ and $\langle m_3\rangle$ with the other integrals. Preliminary numerical experiments confirm the practical usability of these formulas in order to recover the moment of magnetizations. A research report is currently being written to sum up these results.

Finally, a third more ambitious approach has been investigated. As an attempt to generalize the previous expansions, our initial question was: given measurement of $B_3$, say on a square, find a function $\phi (x,y)$ such that $\int \int \phi (x,y)B_3(x,y)\mathrm{d}x\mathrm{d}y$ is a the best possible estimate of the net moment components $\langle m_i\rangle$ ($i=1,2,3$). This problem does not admit a solution because, for any $ϵ>0$, there exists a function ${\phi }_{ϵ}$ allowing to estimate the moment with an error bounded by $ϵ$. However, when $ϵ$ tends to zero, the function ${\phi }_{ϵ}$ is expected to have strong oscillations, which hinders an accurate computation of $\int \int \phi (x,y)B_3(x,y)\mathrm{d}x\mathrm{d}y$ since $B_3$ is only known on a discrete grid of points. We therefore expressed the problem as a bounded extremal problem (see Section 3.3.1 ): to find the best ${\phi }_{ϵ}$ (with the smallest possible error value $ϵ$) under the constraint that $\parallel \nabla {\phi }_{ϵ}{\parallel }_2\le M$. Here, $M$ is a user-defined parameter. We proved theoretical results regarding this bounded extremal problem (existence and uniqueness of a solution, characterization of its solution as a solution of integro-differential equation) and we are currently designing a numerical procedure to compute it. An article on this topic is in preparation.

Still in the course of D. Ponomarev's PhD research, the study of a 2D spectral problem for the truncated Poisson operator in planar geometry has been pursued. It is a simplified formulation of the relation between the magnetization and the magnetic potential (of which the magnetic field is the gradient) and is expected to produce an efficient representation basis (the eigenfunctions of the magnetization-to-field operator). This is a long-standing problem. Noteworthy properties of solutions have been obtained through connections with other spectral problems and asymptotic reductions for large and small values of the main parameters (distance $h$ from the measurement plane to the sample support and sample support size), yielding approximate solutions by means simpler integral equations and ODEs.

The year 2015 was the last of our “équipe associée” Impinge with the MIT and Vanderbilt University. The final report is available on the web page of the associate team (http://www-sop.inria.fr/apics/IMPINGE/ ). This collaboration is currently supported in part by a MIT-France seed funding from the US side, and we applied for a three-years extension of the associate team.

Inverse magnetization issues from sparse spherical data

The team APICS is a partner of the ANR project MagLune concerning Lunar magnetism, associated to the Geophysics and Planetology Department of Cerege, CNRS, Aix-en-Provence (see Section 9.2.2 ). Measurements of the remanent magnetic field of the Moon let geoscientists think that the Moon used to have a magnetic dynamo for some time, but the exact process that triggered and fed this dynamo is not yet understood, much less why it stopped. In particular, the Moon is too small to have a convecting dynamo like the Earth has. The overall goal of the project is to devise models to explain how this dynamo phenomenon was possible on the Moon.

To this end, the geophysicists from Cerege will go to NASA to perform some measurements on samples brought back from the Moon by Apollo missions. The samples are kept inside bags with a protective atmosphere, and geophysicists are not allowed to open the bags, nor to take out the samples from NASA facilities. Therefore, measurements must be performed with some rudimentary instrument and our colleagues from Cerege designed a specific magnetometer. This device allows them to obtain measurements of the components of the magnetic field produced by the sample, at some discrete set of points located on disks belonging to three cylinders (see Figure 4 ).

This collaboration started this year and some preparatory work was necessary fix conventions used by our colleagues from Cerege in order to handle their measurements. During his Master 2 internship, Konstantinos Mavreas has developed a method based on rational approximation, using the same ideas as those underlying the FindSources3D tool (see Sections  3.4.2 and  7.1.1 ), for the case where the field produced by the sample can be well explained by a single magnetic dipole, whose position and moment are unknown. See his report (http://www-sop.inria.fr/members/Konstantinos.Mavreas/main.pdf ). Konstantinos Mavreas is now engaged in a PhD within APICS and will extend these results to the case of several dipoles.